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An important class of approaches to the description of electronic transport through molecules and 
quantum dots is based on the master equation. We discuss various formalisms for deriving a master 
equation and their interrelations. It is shown that the master equations derived by Wangsness, 
Bloch, and Redfield and by Konig et al. axe equivalent. The roles of the large-reservoir and Markov 
approximations are clarified. The Markov approximation is traced back to nonzero bias voltage and 
temperature, whereas interactions and the corresponding rapid relaxation in the leads are shown to 
be irrelevant for the transport under certain conditions. It is explained why the T-matrix formalism 
gives incomplete results except for diagonal density operators and to second order in the tunneling 
amplitudes. The time-convolutionless master equation is adapted to tunneling problems and a 
diagrammatic scheme for generating arbitrary orders in the tunneling amplitudes is developed. 
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I. INTRODUCTION 

Most approaches employed for the description of tun- 
nehng through molecules and quantum dots fall into one 
of two conceptual classes: In the first, one focuses on 
the dynamics of individual electrons tunneling through 
the system. Their dynamics is often described with the 
help of single-particle nonequilibrium Green functionSfii^ 
This approach is sometimes combined with ab-initio 
mcthods"^'^ or perturbation theory for the interactions on 
the dot. In the second, one focuses on the time evolution 
of the many-particle state of the dot and integrates out 
the effect of the leads. This allows one to treat the strong 
interactions on the dot exactly. Typically, the reduced 
density operator of the dot is considered. Its equation of 
motion is the master equation (ME). The two approaches 
become equivalent in the absence of interactions. 

The ME approach involves two steps. First, one de- 
rives the ME from the von Neumann equation for the 
full system. It tells one how the reduced density op- 
erator changes, based on its present and often its past 
values. It is clear that this requires additional assump- 
tions, since the dot is coupled to the leads, which can be 
in any state. It is of course desirable to get by with only 
weak assumptions. 

Second, one solves the ME to find the time evolution 
or the stationary state. Finding the time evolution is 
more complicated if the ME contains the history of the 
dot. It is thus desirable to obtain a ME that is local in 
time. There exist both approximate and exact methods 
for achieving this. 

The ME approach comes in a number of fla- 
vors, among them the original Wangsness-Bloch-Red- 
field (WBR) approachf^i^ii the superoperator formal- 
ism, the T-matrix formalism, the Keldysh-contour for- 
mulation of Konig et al.^-'^i^^^^ and the so-called time- 
convolutionless (TCL) MEii^ The purpose of this paper 
is to clarify the interrelations between these different pic- 
tures and to analyze some of their problems. It is hoped 
that this will facilitate the comparison between results 



obtained with different methods. In addition, the TCL 
formalism is generalized to the tunneling case and is ar- 
gued to provide a powerful tool for studying the dynamics 
of a dot under bias. 



II. THE WANGSNESS-BLOCH-REDFIELD 
MASTER EQUATION 

A. Conventional derivation 

This approacb^ii is commonly described in text- 
booksj^ii^ii^ii^ though not in relation to particle trans- 
port. Several groups have recently applied it to tun- 
neling through molecules We start from 
a Hamiltonian H = Hdot + -ffioads + ^fhyb, where the 
terms describe the dot, the leads, and hybridization be- 
tween them, respectively, -ffdot may contain vibrational 
or spin degrees of freedom and their coupling to the elec- 
trons. The time evolution of the density operator p of 
the full system is described by the von Neumann equa- 
tion p — ~i [H, p], where % = \. We wish to find the ME 
for the reduced density operator pdot{t) = tricadsp(i), 
where the trace is over many-particle states of the leads. 

The central assumption is that i/hyb can be treated 
perturbatively. Operators A are transformed into the 
interaction picture with respect to i/hyb, 
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The density operator in the interaction picture satisfies 
the equation 

p/ = -i [Hhyb,/,P/]- (2) 

Integrating this equation from <o to t yields 

Pi{t) = P/(io) -i I dt' [Hi,y^j{t'),pi{t')]. (3) 

Jta 

Inserting this again into Eq. ^ one obtains 
pj{t) = -i [i/hyb,/(t),P/(io)] 
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dt' [i/hyb./W,[i?hyb,/(i'),P/(i')]]- (4) 

to 

Continuing the iteration, one generates equations con- 
taining arbitrary powers of i/hyb- However, in all of them 
all terms except for the one with the highest power con- 
tain Pi only at time tg. This is problematic when we 
introduce approximations. For example, if we describe 
cotunneling (fourth order), sequential tunneling (second 
order) would only appear in the dependence on the initial 
conditions. Section Hi Bl explains how to avoid this. 

We now consider the initial condition that the system 
was in a product state at an early time Iq, 



pita) = Pdot(io) ^pLds 



(5) 



with PiJjads describing separate thermal equilibria of the 
two leads; their chemical potentials and temperatures can 
be different. A product state is equivalent to dot and 
leads being statistically independent at time ^Oi which 
is natural if iJhyb is switched on at time to. The initial 
condition ^ implies that 



tricads^^hyb,/(i) Pl{to) — 0, 



(6) 



since i?hyb is a sum of terms containing a single lead- 
electron operator. Thus each term changes the lead elec- 
tron number and thus gives zero under the trace, since 
pJoads oi^ly contains states with sharp electron number. 
The trace over lead states of Ec^. (|3|) is 



Pdotj(t) - - / di'trieads [Hhyh.lit), [Hi,ykj{t'),pi{t')]] . 

(7) 

The first term in Eq. ^ drops out due to Eq. ©. Up 
to this point, the results are exact. 

The integral in every term in Eq. ([7]) is of the form 



trioads 

n Pleads,//?/) n' 



dot,/, 



(8) 



where the first and last factors are products of zero or 
more dot electron operators and Y\ i?icads,/ is a product 
of two lead electron operators. The operators may have 
different time arguments. At this point it is usually ar- 
gued that the tunneling should have negligible effect on 
the leads, since they form a large reservoir?^ Therefore, 
in any term one makes the replacement^ii^ii^ 

trieads W B\cads,lPl ~ trioads W -BlcadsjPdot,/ ® Pleads 
= Pdot,/ ® trioads W Pleads jP^oads- (9) 

Here, one replaces any two-time correlation function of 
the leads by the correlation function in equilibrium. Gar- 
diner and ZoUer— point out that one only has to make 
this assumption in the second-order terms. 

In fact, we must only make it in the second-order 
terms: If we were to argue that since tunneling should 
have negligible effect on the leads we can replace p{t) by 



Pdot (i) ® pfoads globally in the von Neumann equation, we 
get trivial results. By taking the lead trace, we obtain 

Pdot = -i trioads [-ff, PdotfX" Pleads] = tricads [^^dot , Pdot] 

/"loads' "wtiich is just the unperturbed time evolution. This 
is exact if p{t) — pdot{t) (X" Pleads holds, but this is not 
very useful, since the condition is generally not satisfied 
at any other time. We will see that in the superoperator 
approach we do not have to worry about this, since we 
only assume a product state at an initial time to, as in 
Eq. ([5]). Furthermore, in the TCL approach we can avoid 
even this assumption. 

If approximation ^ holds, Eq. ([7]) becomes 

Pdot,l{t) — dt'tricads [Hhyh.l{t), [Hhyhj{t'), 

Jto 

Pdot,/(t')®pLds]]- (10) 

This ME is non-local in time. 

To make it local, one usually introduces the Markov ap- 
proximation, which replaces Pdot,/(^') by Pdot.i{t)- This 
means that the rate of change of pdot,/ at time t is de- 
termined by pdot,7 at the same time t only. The approx- 
imation is usually motivated by an argument of the fol- 
lowing typeji^ii^ Eq. PU)) contains two-time correlation 
functions for the leads of the form ^ . These correlation 
functions decay rapidly on the time scale of the dot dy- 
namics so that they can be replaced by (5-functions. We 
come back to this point in Sec. Ill Bl 

The same assumption also implies that as long as t — to 
is large compared to the lead correlation time, one can 
replace to by — oo. With t' = t — T one obtains 

Pdot.l{t) = - rfr trioads [i?hyb./(0 J^hyb./(i - t), 

Jo 

Pdot,7(0®PLds]]- (11) 

Transforming back into the Schrodinger picture using 
Pdot,7(t) = e*^<i°'*pdot(t)e-'-f^<'°** one finds 



Pdotit) ^ -i[Hdot, Pdot{t)] - / drtrioads 

Jo 



hyb. 
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Pdot(t) f^pLds 



(12) 



which is local in time. The first term describes the unper- 
turbed time evolution of pdot and the second is a correc- 
tion of second order in iJhyb- The restriction to second 
order entered when we made the large-reservoir approx- 
imation after iterating the equation of motion to second 
order. More explicit expressions are given in App. [X] 
They also show that, if the tunneling amplitudes do not 
depend on wave vector, the leads only enter through their 
density of states, temperature, and chemical potentials, 
regardless of interactions in the leads. 

At this point it is often assumed that the off-diagonal 
components of pdot decay rapidly and can be neglected. 
For some components this can be motivated by super- 
selection rulesiS^i^"^ If two dot states |m), \n) differ 
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in an observable that couples strongly to the environ- 
ment, unavoidable interactions lead to rapid decay of 
superpositions of these states and thus of /o^n-" The 
standard example is the charge.-—'^-^ Due to Gauss' law, 
the effect of the charge can in principle be measured 
equally well on any arbitarily large surface surrounding 
the systemi^ Therefore, superpositions of dot states with 
different charge are not observed. 

On the other hand, the description of spin precession^^ 
requires the off-diagonal components. Different spin 
states also differ in their long-range (dipole) fields, but 
these fall off more rapidly than the Coulomb field. This 
suggests that Gauss' law^l is crucial for superselection 
rules and not just any algebraic decay. 

If all off-diagonal components decay rapidly, one is left 
with the diagonal components Pm = Pmm' i-^-i prob- 
abilities of dot many-particle states \m). The principal- 
value terms in Eq. (|A1[) then cancel and one obtains 



Prr 



-2^EEl«*l("''l^hyb|p)b-»r 

ij p 

X {W, P,n - W, Pp) 5{Ep + e, - E„ 



(13) 



Here, Wi = ((« | Pleads I *)) ^^'^ probability to find the leads 
in state |z)). Defining the transition rates 



^T^,|((z|(m|77hyb|n)b-))| 

X S{Ep + ej - E„i -Ei), 
we obtain the well-known rate equations 



Prr 



(14) 



(15) 



The first term describes transitions from other states |n) 
to state |m), whereas the second describes transition out 
of state I to). The rate equations imply local conservation 
of probability — Pm only changes due to probability flow- 
ing into or out of state |to). This conservation law can 
be implemented in a gauge theoryi^ 



B. Discussion of the Markov approximation 

The Markov approximation is usually motivated by 
rapid decay of the lead correlation functionsJ^'i^'^ In 
the second-order approximation, each nonvanishing term 
contains one creation operator and one annihilation 
operator a. The result is non-zero only if both belong 
to the same single-particle state. The trace tricads over 
lead many-particle states is replaced by a sum over single- 
particle states characterized by lead index a, wave vector 
k, and spin a. As discussed in App. \^ the correlation 
functions are Green functions , G-* . If the leads are 
normal metals, these decay on the time scale of the quasi- 
particle lifetime. (The non-quasiparticle background in 
the spectral function is broader than the quasiparticle 



peak, corresponding to faster processes, which are less 
critical for the validity of the Markov approximation.) 
However, the quasiparticle lifetime becomes long at the 
low temperatures at which experiments are performed. 
Does the Markov approximation break down in the ex- 
perimental temperature range? This question is also rel- 
evant because the lead Hamiltonian -ffioads used in actual 
calculations does not contain any interactions. Thus in 
our model, there is not broadening of the quasiparticle 
peak and , G^ do not decay in time. 

While each term separately does not decay, their sum 
does. We replace the sum over k by an integral over 
energy, including the density of states. At low tempera- 
tures we can restrict the integral to the energy window 
between the two chemical potentials /i<, /i>. Assuming 
a constant density of states and k-independent tunneling 
amplitudes, we end up with integrals of the type 



dEe 



±iET 



±- 



(16) 



This expression contains a typical time scale h/eV = 
Tieads for the decay of correlations, restoring Planck's 
constant for the moment. Thus the energy governing 
the decay of correlations is the bias, not the electron- 
electron interaction. The same energy scale determines 
dephasing, i.e., the decay of superpositions due to differ- 
ent chemical potentials in the leads. For arbitary tem- 
peratures, the limits of integration are roughly /i< — fc^T 
and fiy+ ksT and the characteristic time is the smaller 
of h/eV and h/kBT. Note that the contribution from 
the quasiparticle lifetime is proportional to^ 1/T^ and 
is thus irrelevant at low temperatures. 

For weak tunneling, or specifically if the conductance 
is small compared to the quantum conductance. 



(17) 



the typical time between two tunneling events is tq = 
e/I ^ h/eV = Tieads- Then the dot dynamics is indeed 
much slower than the decay of lead correlations and the 
Markov approximation is justified. It follows from the 
weak-tunneling approximation, which we have to make 
in any case to work in low-order perturbation theory. 

This argument may fail if tunneling events are strongly 
correlatedi^2i2ii2^i2& In this case, two or more tunneling 
events can often take place during a time much shorter 
than To = e/I and the relation (fT7| does not guarantee 
the validity of the Markov approximation. 

A related point is seen if we proceed slightly differently 
in the derivation. Starting from Eq. ^ and inserting 
Eq. ([3]) with renamed variables 



pi{t') = pi{t) - ^ j dt" [-ffhyb,/(i"), Pi{t'% (18) 



we obtain 



-i [Hhyhj{t),Plito)] 
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- / dt'[i?hyb./W,[^fhyb,/(i'),P/W]] 
Jto 

+ 1 f dt' f dt"[Hi,ybj{t),[Hi,ybj{t'), 

Jto Jt 

[ffhyb,/(i"),P/(i")]]], (19) 

which is still exact. If we now restrict ourselves to the 
second order in i^hyb, we can drop the last term. We 
have obtained an equation that is local in time without 
invoking the Markov approximation. 

We have pushed non-local terms into higher orders in 
-ffhyb- Iterating the procedure, we can achieve this to any 
order. All relevant terms contain pi{t) instead of p/(to), 
which appears in the naive expansion in Sec. Ill Al Thus 
we can for example derive a ME containing sequential 
and cotunneling contributions. 

Now we can make the large-reservoir approximation 
and replace pi{t) by Pdot,i{t)® P\cs,ds^ as above. But here 
we perform this replacement only at time t. We then 
obtain, by tracing over the leads, 

Pdat,l{t) — ^ I dt' iV\cads[H\iyh,l{t),[Hhyh,l{t')i 

Jto 

Pdot,/(i)®pLds]]- (20) 

We finally replace by — oo. This is a remnant of the 
Markov approximation, but only if the "true" is finite. 
The choice of can be viewed as a part of our model 
as opposed to the approximations employed to solve it. 
Equation (|20|) can also be obtained from a variational 
principle, again only expanding in i?hyb without explicit 
Markov assumption)^ 

We have obtained the same local ME in seemingly dif- 
ferent ways. The explanation is that we have not made 
independent approximations. The rapid decay of the lead 
correlation functions follows from the assumption of weak 
tunneling, which also allows us to use perturbation the- 
ory in -ffhyb- On the other hand, the large-reservoir (or 
Born^i) approximation pi{t) = PAot,i{t) ® Pleads logi- 
cally independent. 



III. SUPEROPERATOR FORMALISM 

The WBR ME can also be derived in the superoperator 
formalismfi^i^i^i^S which facilitates expansion to higher 
orders in iJhyb- We here define a superoperator as an 
operator acting on the space of linear operators on the 
Hilbert space. The von Neumann equation is written as 
p = —iCp, where 

^ = ^dot + Pleads + -Chyb (21) 

is the Liouvillian defined by LdotP = [Hdot, p] etc. 
The solution reads pit) = e~*^(*~*°V(*o)- We define 
projection (super-) operators^ML Q by Vp^t) = 
[tricadsP(i)] ® Pleads ^'^'^ Q = ^ - 'P ■ Notc that V maps 



a density operator onto one in product form with the 
leads in equilibrium, while retaining the information on 
the dot state. Conversely, Qp contains the information 
on the leads and the dot-lead correlations. It is easy to 
prove the identities 

VCdot = CdotV, (22) 

'P^loads = •Clcads'P = 0, (23) 

rCi-^ybV = 0. (24) 

The last relation is essentially equivalent to Eq. ([5]). 

The von Neumann equation can be split into two 
equations^ii^i^i^a 

^■Pp = -iVCVp-iVCQp, (25) 
^Qp - -iQCQp-iQCVp, (26) 
which can be solved by Laplace transformation, 

/•oo 

F(s)= / dte-''F{t), (27) 

JQ 

where we have set t^ ~ Q. We find 

sVp'Vp{0) = -iVCVp-iVCQp, (28) 
sQp-Qp{0) ^ -iQCQp-iQCVp. (29) 

Inserting the solution of the second equation, 

Qp = {s + iQCy'^QpiO) -i{s + iQCy^QCVp, (30) 

into the first, we obtain 

sVp = Vp{Q)-iVC{s + iQC)-^Qp{{))-iVCVp 

-VC{s + iQC)-^QCVp. (31) 

Using Eqs. JUll-dlll), we find 

sVp^Vp{Q) 

- iVCuyh [s + iCdot + iCicads + iQChyhT^QpiO) 

- iVCdotVp (32) 

- VCiiy], (s + iCdot + «Aoads + «QAiyb)~^QAiybPP- 

We can insert another Q following Q£hyb since Q = 
QQ. This makes the resulting expression £dot + £ioads + 
QChyhQ hermitian. 

Transforming back into the time domain and shifting 
the initial time back to to we find 

^Vp = -^7'/:hybe-*(^^°'+^'-^=+S^''-•=2)(t-to)gp(^g) 
-iVCdotVp{t)-VCY,y], f dt' 

Jto 

X g-»(£dot+£lcads + Q/:hybS)(t-t')£j^yj^7?p(t'). (33) 

(We will show later that the projections Q in the expo- 
nentials remove all reducible terms from the expansion 
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in powers of £hyb-) Equation (j33p is an exact ME, which 
is non-local in time. 

Starting from Eq. ((33|) , the weak-coupling limit as dis- 
cussed in Ref. [l^ now consists of (a) neglecting all pow- 
ers of £hyb beyond the second and (b) dropping the de- 
pendence on the initial condition for Qp(io)- (a) is just 
the weak-tunneling approximation of Sec. Ill Al (b) ne- 
glects a term of linear order in £hyb and has to be shown 
to be consistent. The rationale given in Ref. [H is two- 
fold: First, the term in Qp{tQ) is a small (of order iJhyb) 
correction to Vp^to), and second, it is not accumulated 
over time, being a correction to the initial conditions. 
These arguments appear to be weak: Qp{to) and Vp{to) 
lie in orthogonal subspaces and it is not obvious that 
their magnitudes can be meaningfully compared. Fur- 
thermore, Eq. (|33p shows that Qp{to) does affect 'Pp{t) 
for all t > to, even to first order. 

Dropping the dependence on Qp{to) is trivial if 
Qp{to) — 0. This is not an approximation but an ini- 
tial condition, see Sec. Ill Al 

With approximations (a) and (b), Eq. ([33]) becomes 



dt 



Vp - ^iCdotVpit) 



■P/Ihyb / di'e^*(^''°'+^'-<'=)(*-*')£hyb'Pp(t')- (34) 

Jtn 



Inserting the definition of V and writing the Liouvillians 
as commutators we obtain 



^ Pdot ^ -i[i?dot, Pdot(i)] 

- / di'trieads[ifhyb,e-''(^^-+^'-^=)(*-*')[i?hyb, 



Pdot(i')'»p"cads] 



(35) 



Now for any Hamiltonian H with associated Liouvillian 
C and any operator (not superoperator) A the identity 



(36) 



holds. In our case we thus get 



-77 Pdot — —i[Hdot, Pdot{t)] — / C?i' tricads[-f^hyb, 

rg-i(-ffdot+Hload=)(t-t')^ gi(-ydot+ffload=)(t-t') 



hyb < 

e-^^-'(*-*'Vdot(t')e^^-'(*-*')^pLds] 



(37) 



Compare this to the WBR result, Eq. To get 

there, we have to replace e~'^<i°'(*~* ■'pdot(i')e*^''°'^*~* ^ 
by Pdot{i)- This is nearly the same: pdot{t) is described 
by the full Hamiltonian H , while Eq. ([37|) only contains 
the unperturbed time evolution due to -ffdot . This is con- 
sistent with the second-order approximation, since any 
correction to the unperturbed time evolution adds more 
powers of i?hyb- Thus the Markov property again follows. 
We can extend the range of integration to to ^ ~oo ar- 
guing as in Sec. IIIBI 



Note that we did not need the large-reservoir assump- 
tion p{t) = Pdat{i) ® Pleads ^^^^ approach. We have 
only assumed the density operator to be of this product 
form at an early time — a much weaker assumption. 



IV. THE T-MATRIX APPROACH AND 
FERMI'S GOLDEN RULE 

The T-matrix approacli22^'25>3642^44 
its leading-order approximation, Fermi's Golden 
Rule^^i^i^ are used by several groups to describe 
tunneling processes, since they provide a straightforward 
derivation of the transition rates in the diagonal rate 
equations. The approach is presented in many text- 
books. Bruus and Flensberg^^i discuss it in relation to 
the tunneling problem. 

The derivation starts out by writing 



H[t) = ffdot + H, 



leads 




(38) 



where ry is small and positive. Thus the hybridization is 
switched on very slowly. We assume that the system was 
in an eigenstate \i) of Hq at time to- The probability that 
it is in another eigenstate |/) at time t reads 



mm? 



(/ITexp 



rt N 2 

-I j dt'Vi{t') I 12) 



(39) 



where T is the time-ordering operator. The transition 
rate between states \i) and |/) is then defined as 

r;,^||(/|z(t))p. (40) 

Taking the limit rj O'^ and defining the T-matrix 

1 



hyb ■ 



hyb 



E, -Hn+ iO+ 



1 



hyb 



Ei — Hn 



iO+ 



1 



hyb 



E,-Ho + iO+ 



hyb 



one obtains 



Tj, = 2n5{E.,~Ef)\{f\T\i) 



(41) 



(42) 



where Ei and Ef are the eigenenergies of states \i) and 
I/), respectively. The leading order is obtained by re- 
placing T by iJhyb, 



r/, = 2^<5(i?,~i?;)|(/|i/hybN) 



(43) 



which is Fermi's Golden Rule. 

To draw the connection with the rate equations, we 
choose the states \f) as product states of many- 
particle states of the dot, |m), and of the leads, \ij). 



6 



Summing Tfi over the lead states we obtain the tran- 
sition rates 

Rn^rn = 2^ ^ W^, | ((/ 1 (TO|T|n) | ',5(£;„ + e, - ~ £/) 

(44) 

from state \n) to |m). Here, E^a (q) are eigenenergies 
of dot (lead) states and Wi is the probability to find the 
leads in initial state \i) at time io ~* To write 

down Eq. , one has to make the assumption that the 
probability Wi is independent of the state of the dot. 
This means that the system is in a product state p = 
Pdot ® Pleads at time t^. This is the same assumption 
usually made in density-operator approaches. If the leads 
are in equilibrium at time Iq one can express Wi in terms 
of Fermi functions. 

In the next step, i?„_>„i is identified with the transition 
rate i?„^„j appearing in the rate equations. However, 
what is actually calculated is the rate of change of the 
probability of state |to) under the condition that the dot 
was in state |n) at time io — oo, cf. Eqs. (|40l) and (|44| . 
On the other hand, in the density-operator approach one 
calculates the rate of change of the probability of state 
|m) under the condition that it is in state \n) at the same 
time immediately before a possible transition. The two 
are the same only if the dot remains in state \n) from time 
to through t. This is of course not usually the case. 

In the sequential-tunneling approximation, one evalu- 
ates the rates to second order in -ffhyb- Since two powers 
of 77hyb are required for the final transition to state |to), 
no transitions can occur between io and t. Therefore, to 
second order, where the T-matrix approach reduces to 
Fermi's Golden Rule, it gives the same transition rates 
Rn^m as the WBR approach. Beyond leading order, 
Rn^m and Rn^m dcscribc different quantities. 

To leading order, T ~ ^^hyb, the rates Rn-,m m 
Eq. are indeed identical to the the WBR result, 

Eq. (fTil) . The latter has been obtained under the Markov 
assumption. One might wonder where the Markov as- 
sumption entered in the T-matrix formalism. It is im- 
plied in the derivation, since to second order p does not 
change between times and t anyway. 

The rate equations ([TS]) appear to be obvious. How- 
ever, can they be derived in the T-matrix framework? 
Certainly not without further assumptions, since they 
omit the off-diagonal components of pdot necessary for a 
complete description. 



V. THE KELDYSH-CONTOUR APPROACH OF 
KONIG ET AL. 

Konig et aZ4i2i i^°i^^ have developed a diagrammatic 
technique to generate a perturbative expansion in the 
tunneling amplitudes. This approach has also been ap- 
plied to tunneling through molecules j^^'^^ In the present 
section we show how it is related to WBR theory. Ref- 
erence 1^ concerns a quantum dot with electron-electron 



and electron- vibration interactions. A unitary transfor- 
mation replaces the latter with an exponential operator 
in the tunneling Hamiltonian4^ We do not consider this 
transformation here. This does not restrict the models 
covered by our discussion — we include any bosonic modes 
and electron-boson interactions into Hdot- 

As usual, the system is assumed to be in a product 
state at an early time to with the leads in separate equi- 
libria. The propagator H of pdot from time to to t > to 
is defined by 

Pdot(t) = H(t, to) pdot(to). (45) 

(We use the same notation and time ordering as elsewhere 
in this paper.) H(t, to) is represented by a diagram on the 
Keldysh contour between times to and t. Fig. 3 in Ref. 
Konig et al. then identify the irreducible part of H (a 
subscript is added to distinguish T,k from E introduced 
above), defined as the sum of all diagrams that cannot be 
cut at an intermediate time without cutting a lead line 
representing the pairing of ajj^.^. and aako- • H is expressed 
in terms of E/f by a Dyson-type equation, 

H(t,to) - H(o)(t,to)+ / dt2 f ' dhn^"\t,t2) 

X EK(t2,ti)H(ti,to), (46) 
containing the bare propagator 

[Tl^'Ht^to)]^^, = <5w<^™'e-^(^"-^")(*-*''). (47) 

In Ref. 9, Eqs. (|l5|) - ([T7|) are presented for arbitary initial 
time t' > to. Equation ([45]) would then imply H(t,t") = 
H(t, t') H(t', t"), which contradicts Eq. This is not 

a problem since the equations with initial time to are 
sufficient for the derivation of the ME. 

We want to find the propagator H explicitly. Clearly, 
we have Vp{t) — p^-'-C-i't^^o) p{tf^'). Since the initial con- 
dition Qp{to) = is assumed, we find 

Vp{t) = [T'e-'^^*-*") V] Vp{h) = H(t, to) Pp(to). (48) 

This defines the propagator H for initial time to. The 
last factor in H is expendable. 

Can we find Y^k to satisfy Eq. We expand the 

exponential in Eq. (pS)) , noting that all lead creation and 
annihilation operators must be paired for the result to 
be nonzero since P£hyh'P = 0. Diagrammatically, this is 
represented by a lead-fermion line connecting two i?hyb 
insertions 1^ The lead trace then gives a Fermi factor for 
any pair of insertions. With regard to this pairing, we 
can identify the irreducible part S^f and write 

Pg-*c(t-to)-p ^ -pg-«(z:dot+Cioad=)(t-to) / / dti 

Jto Jto 

X p-i(£do,+£lcads)(t-t2) X;^(t2,tl) 
X g — *(-Cdot+£^lcada)(tl— to)-p _j_ 

^ -Pg-i(£do,+£lcads)(t-to) + 7? /" dt2 dtl 

Jto Jto 

X g-»(£dot+£ieads)(t-t2) 5]^(^2,ti)e-*^(*i-*°)7'. (49) 
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If we insert V at some intermediate time tn in Eq. ([491) 
this forces all lead-fermion operators to be paired for 
times smaller than and larger than tn separately. This 
is equivalent to the diagram being reducible at time i„. 
Conversely, if a diagram is reducible at time tn , inserting 
V there does not change the result. Consequently, 

J to J to 

X e-^(^i°'+^i-''=)(*-*2)SK(i2,ii)^e^*^(*i-*o)7'. (50) 

Here, we can identify the bare propagator 

n(°)(t,t') = 7?e~*(^<'-+^'-''=)(*-*'). (51) 

These two equations correspond to Eqs. (|46p. (I47p except 
that we had to introduce another projection V at the 
final time. Apart from this, we recover the Dyson-type 
equation of Ref. d for initial time to. 

Inserting the Dyson-type equation (jSO]) for the propa- 
gator into Eq. (|48l) and taking the time derivative, one 
obtains a ME^ 

-Vp{t) = -iCdotVp{t) + v / dt'^KMvp{t'). 

dt Jto 

(52) 

This corresponds to Eq. (25) in Ref. 9 (a factor of V on 
the left can be included into the definition of Y.k)- 



Q 



FIG. 1: Diagram for the second term in the ME (|53|l . 
Time is increasing to the right. The bare propagator 
g-»(-Cdot+£ioads)(«-« ) fQj. ^jjg density operator (not occuring 
here) would be represented by just two lines. The full prop- 
agator e 



-iC(t-t') 



part e 



-iS£S(t-t') 



is shown as a hatched box. Its irreducible 
is distinguished by writing "Q" beside it. 
All additional projection operators are likewise indicated by 
"P" or "Q". An insertion of —i£hyh is denoted by a heavy 
bar connecting the two lines. The projected density operator 
Vpit') at time t' is shown as a filled semicircle. 

We can gain further insight by returning to Eq. 
restricted to the case Qp(io) — 0, 

d /■* 

— Vp^-iUotVp{t)~VC^y], / dt' 
at Jto 

Comparing to Eq. ([5^ we find 

(54) 

Thus the ME of Ref. is equivalent to the WBR ME to 
all orders in -ffhyb and we have derived an explicit ex- 
pression for the irreducible part. For later, we introduce 
in Fig. [1] a diagrammatic representation for the second 
term in Eq. (|53p . 



Now the factors of Q in the exponential find a natural 
interpretation: They remove all reducible terms from the 
expansion of Eq. ((M)) . This is because at a point where 
a diagram is reducible, one can insert a V. But if a Q is 
present at this point, we obtain VQ — QV = 0. Thus all 
diagrams in the expansion that are reducible to the left 
or to the right of an insertion of £hyb vanish. 



VI. THE TIME-CONVOLUTIONLESS MASTER 
EQUATION 

This approach leads to a ME that is local in time 
and exact, and thus avoids the Markov assumption. The 
Markov assumption is valid for weak tunneling and not 
strongly correlated tunneling events, as discussed above, 
but becomes increasingly dubious at higher orders in 
^^hyb or in resummation or non-perturbative schemes. 
The approach was developed by Tokuyama and Mori^^ 
and others^i^i^ and is discussed in Ref.[l^. This section 
adapts it to the tunneling problem. 

We again start from Eqs. pS)) and (pS)) and express 
Qp{t) in the first equation in terms of Vpit) and Qp{to) 
with the help of the second. We do not make any as- 
sumptions on Qp{to). The solution of Eq. reads 

Qp{t) = e-^S£Q(t-to)g^(i^) 



Jto 



(55) 



Next, we express p{t') by propagating the full density 
operator backward in time, 



Qpit) 



^ g-iQCQ(t-to) 



'to 



'Qp(to) 

X Q£Pe-^^(*'-*) [Vp{t) + Qpit)] . (56) 
Moving all terms in Qp{t) to the left we obtain 

(1 - S) Qp{t) = e-«2^S(*-*o) Qp(t^) + s T?p^t) (57) 
with the superoperator 



S(i-io) 



Jtn 



The time argument t — to will be suppressed if confusion 
is unlikely. The integral can be performed, giving 



E(f - to) = f dt' f e-'S^S(*-*') Q A 
Jto \ dt' 

^_ ^-iQCQ(t-t') Q^-iC(t'-t)^ 

-iQCQ(t-to) Q^-iC{to-t) _ (gg) 

is more suitable for the expansion 



I to 

+ 

Q 



The integral form 
in £hyb, though. 
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Applying the inverse (1 — S) ^ to Eq. ([57)) . we obtain 



(63) 



Qp{t) = (l-S)-ie-'2^2(t-to)g^(i^) + (l_j^)-i5^p^(i), 

(60) 

This remarkable equation asserts that we can reconstruct 
Qp and thus p = Vp + Qp at time t from Vp at time t 
and Q/9 at some arbitrarily early time to, even though 
Vp only contains information on the dot state. 

Inserting Eq. ((60|) into Eq. ([25|l we obtain an equation 
of motion for T'p alone, 



-ivc{i-j:y^Vp{t) 

^ iVC{l - E)-ie-*2£2(t-to)g^(i^)^ (61) 



which, together with Eq. ([55)1 or (|5^ . constitutes the 
TCL ME. It indeed only contains Pp at the time t. 
Eqs. (pTj) and ([55]) can be generalized for time-dependent 
Hamiltonians^i^ by replacing the multiplication with 
t — to by a time integral. 

The TCL ME is exact, but relies on the condition that 
the inverse of 1 — E exists. For a different system not 
involving transport, one can find cases when 1 — E is 
singular!^ Since S vanishes att = to, cf. Eq. ((59)) . Breuer 
et alJ^'^ conclude that 1 — E is not singular for suffi- 
ciently small t — to- This does not directly apply to our 
case, since the time dependence is governed by the dy- 
namics of the full system, which permits arbitarily large 
excitation energies. In other words, the eigenvalues of C 
and QCQ in Eq. ([55)1 are not bounded. 

Using Eqs. (PT|) - ([M)) we obtain a more explicit form. 



dt 



Vp{t) = -iCdotVpit) ~ iVUyb (1 - E)-ipp(t) 

— iVChyh (1 — E)^"*" e~*(-^<iot-|-£leadB-|-Q£hybQ)(*-to) 

X Qp{to) (62) 



with 



We have used that 7^(1 — E)^^ — V, as can be seen 
by expanding in powers of E. Equation (j62p should be 
compared to the non-local ME ([55]) . 

For the case Qp{to) = 0, the equation simplifies to 



- Vp{t) = -iCdotVpit) - iP£hyb[l - E(t - toT^Vpit). 

(64) 

It is then tempting to take the limit to — > —oo. One has 
to check whether this limit exists for E(t — tg). 



J to 



As advertized, the TCL ME ((521) is local in time, al- 
though the dynamics is generally not Markovian. The 
memory effects have been shifted into the time depen- 
dence of the coefficients of 77^ ^55,56,57 ^j^jg works because 
the integro-differential WBR ME is linear. One can then 
show that a purely differential equation with the same 
solution exists. Maniscalco et alP use the exact so- 
lution for the damped harmonic oscillator to illustrate 
that non- Markovian dynamics is indeed compatible with 
a TCL formulation. This and related result s^^'^^ do not 
involve transport. 

We briefly comment on the question of positivity of 
the reduced density operator, i.e., the requirement that 
all its eigenvalues are non-negative. As in Ref. H^, the 
coefficients in the TCL ME are time-dependent. Thus 
Lindblad's-^ criterion for positivity does not apply. How- 
ever, the TCL ME is exact so that its solution for Vp sat- 
isfies Vp ~ \^T:\ca.dsp\® Pleads times and thus certainly 
satisfies positivity. It is an important question whether 
perturbative approximations destroy this property. 

To obtain the sequential-tunneling approximation to 
Eq. (j62p . we expand to second order in £hybj 



dt 



Vp{t) ^ -iCiotVp{t) - T'Aiyb / dt' e-^(^^°'+^i-d=)(*-*')/:hybe"*^'^"°'+^'"-"=^(*'"*^Pp(i) 



7'/:hybe"'^^'^°'+^'-'''=^(*"*''^Q/9(io) -^-Chyb / di'e-*(^'^°'+^'-''=)^*-*')£hybe"*^'^'^°'+^'-'''=^^*'"*'') Qpih)- (65) 

J to 



The first term describes the unperturbed time evolution. 
The two inhomogeneous terms describe the effect of a 
deviation of the state at time to from a product state with 
leads in equilibrium. The third term is the only one of 
first order in the tunneling amplitudes, thus for Qp{to) — 
there are no first-order terms. For readers familiar with 



optical response theory this may seem surprising. We 
briefly discuss first-order terms in App. |Bj 
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A. Perturbative expansion 

Since expansions beyond second order are clearly cum- 
bersome to write down, a diagrammatic representation is 
helpful. We here assume the simplifying initial condition 
Qp{to) = 0. We first expand Eq. ([64)) in powers of E, 



dt 



Vp{t) = -iCdotVpit) 

-i7'/:hyb(E + ES + EEE + ...)-PpW. (66) 



The first-order term vanishes since VChyhV — 0. Note 
that the series 1 + S + -I- . . . might not converge even 
if the inverse of 1 — S exists. 



Lit -Q = Ht 



FIG. 2: Diagram for the superoperator S(t — to). The in- 
terpretation of symbols is given in the caption of Fig. [T] The 
right-most superoperator in Eq. (|63|l corresponds to the lower 
right corner of the diagram. 




(a) 




Q QP 

^^^^^^ 




(b) 



FIG. 3; Contributions to the TCL ME (|66l containing (a) 
one and (b) two powers of the superoperator E. The filled 
semicircle denotes the density operator Vp(t). In (b), no time 
ordering of t' and t" is implied. 



The superoperator S, Eq. (f63|) . first propagates the 
density operator backward in time, projects it, inserts a 
perturbation vChyb, and then propagates it forward again. 
Its diagrammatic representation is shown in Fig. [2l The 
contributions to Eq. (j66p with one and two powers of E 
are represented by the diagrams in Fig. [3) It is obvious 
how the series continues. 




FIG. 4: General form of all terms involving tunneling in the 
TCL ME (|66p . Here, the filled semicircle denotes the pro- 
jected density operator Vp at time t' . 



This is a good place to compare to the approach of 
WBR and Konig et al. The contribution from the tun- 
neling is in this case given by the last term in Eq. (j53p 
or the diagram Fig. [1] To show its consistency with the 
TCL equation, we note that all tunneling contributions 
in the TCL approach are of the form shown in Fig. U) 
This is indeed just a deformation of Fig. [1] (the "Q" ad- 
jacent to the final £hyb is expendable). The upper part of 
S is thus identical to the irreducible part Y,k without the 
final £hyb- This is also seen by comparing the algebraic 
expressions (|54)) and (|63)). 

In order to expand the ME in powers of £hyb, we next 
expand the exponentials as 



X e 



-i(£dot+£lcad3)(tl-t') 
(•t2 



X e-*(^dot+£lcad3)(t2-tl)£ g 



X e 



-i(Cdot+£loadB)(tl-t') 



(67) 



and 



g — i(£dot+£loada + £hyb)(*'-*) _ g — « ( ^dot +£loadB ) ( * ' — * ) 



hyb 



X g — i(£dot+'Cloads)(il— *) _|_ 



(68) 



The odd terms in Eq. ((68)) obtain an additional minus 
sign, since the reversed time order gives an additional 
minus sign for each integral. 

As discussed earlier, any insertion of V forces an ex- 
pression to be reducible at that point and any insertion 
of Q means that the diagram must not be reducible at 
that point. V and Q thus govern the construction of 
non- vanishing diagrams but do not affect their values. 
For the case of a single factor of E the only diagrams up 
to fourth order in >Chyb are the ones shown in Figs. |S{a)- 
(d). The only fourth-order contribution from the term 
with two E is shown in Fig. |5{e) . The second insertion of 
£hyb counting from the lower right comes from the upper 
factor of E due to the right-most "Q" in Fig.|3I^b). 

The diagrams in Figs. O^d) and (e) are topologically 
equivalent but differ in the time ordering. Diagram (d) 
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FIG. 5: Diagrams for all terms in the TCL ME up to fourth 
order in the tunneling amplitudes. The dashed lines denote 
the pairing of lead-electron operators a, with the same 
quantum numbers. 



has t = t4 > ti > t2 > t^, whereas (e) has i4 > ti and 
ti > t2 > t^. In addition, (e) has an odd number of 
'Chyb on the lower branch of the upper E, leading to an 
additional minus sign from Eq. ([68]) . If we add (d) and 
(e), contributions with ti > t2 cancel and we obtain ^4 > 
t2 > ti and t2 > ta and an overall minus sign. Higher- 
order Diagrams are constructed in the same manner. 



VII. SUMMARY 

In this paper various approaches to the ME for tun- 
neling through molecules and quantum dots have been 
discussed and compared. The standard derivation of 
the WBR ME relies on two assumptions: (1) Weak 
tunneling — allows one to use low-order perturbation the- 
ory in the tunneling amplitude and generally implies the 



Markov assumption of temporal locality. (2) The leads 
form large energy and particle reservoirs — together with 
the first assumption and an initial density operator at 
time of product form with the leads in equilibrium this 
allows one to write the density operator as a product at 
all times. Assumption (2) must only be made when cal- 
culating expressions of the desired order in -ffhyb- Making 
it globally leads to trivial dynamics. 

The Markov approximation does not require one to in- 
voke rapid relaxation due to electron-electron or electron- 
phonon interaction in the leads. In fact, under rather 
weak conditions, these interactions are shown to be ir- 
relevant for the tunneling. A short time scale describing 
the decay of correlation in the leads emerges naturally, 
given by the inverse bias or the inverse thermal energy, 
whichever is smaller. At low temperatures, this domi- 
nates over the contribution from the finite quasiparticle 
lifetime. A more careful analysis is needed if tunneling 
events are highly correlated in time. 

The superoperator derivation of the WBR ME clarifies 
the role of the Markov approximation. No assumptions 
beyond weak tunneling and an initial density operator 
of product form are made. The large-reservoir assump- 
tion (2) is not required, beyond this initial condition. 
The resulting WBR ME is non-local in time. An explicit 
Markov approximation is required to make it local. 

The ME of Konig et al^ is equivalent to the WBR ME 
to all orders in tunneling. Its memory kernel, given as 
a diagrammatic perturbation series in Ref. has been 
written down in superoperator form. 

The T-matrix approach only gives rates for the diag- 
onal components of the reduced density operator. The 
large-reservoir assumption of a product state enters in the 
guise of statistical independence of dot and lead states. 
The T-matrix approach gives the same expression for 
these rates as WBR theory with the Markov assumption 
only to leading (second) order in i?hyb, corresponding to 
Fermi's Golden Rule. This is because in the T-matrix 
approach one calculates a subtly different quantity than 
in the WBR approach. 

The TCL ME^^ has been adapted to the transport 
problem. It is an exact equation for the dynamic reduced 
density operator that is local in time but does not require 
a Markov assumption. It thus works to arbitrary orders 
in perturbation theory. The superoperator E(t— io) play- 
ing a pivotal role in the TCL ME has been given in an 
explicit form. 

A number of technical questions regarding the TCL 
formalism applied to transport are still open. First, un- 
der what conditions does the inverse of the superoperator 
1 — E(t — to) exist, which appears in the ME? Second, 
if we assume the system to be in a product state with 
leads in equilibrium at time t^^ can the limit ^ —00 
be taken? And third, is positivity of probabilities satis- 
fied in perturbative approximations to the TCL ME? Of 
course, this may be asked for any approach. We leave 
these questions for future work. 

The assumption of a product state with leads in equi- 
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librium at time to is not required in this approach, but 
simphfies the results. There are processes of first order 
in iJhyb that are physically reasonable but are omitted if 
we make this assumption. They can be incorporated by 
choosing non-standard initial conditions. 

A diagrammatic scheme for generating arbitrary or- 
ders in iJhyb in the TCL ME has been developed. The 
relation to the diagrams of Konig et al.^ and thus to the 
WBR ME is easily seen. Our diagrams have interesting 
additional structure, since the projected density opera- 
tor is propagated backward in time to make the equation 
local. The diagrammatic expansion to fourth order is 
shown explicitly. 
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APPENDIX A: THE WBR MASTER EQUATION 
TO SECOND ORDER 

Equation ([12]) is the WBR ME to second order in i?hyb- 
We here give two more explicit forms. First, we derive a 
fully general expression useful for later comparisons. We 
introduce dot states |m) and lead states |«)) with eigenen- 
ergies Em and e^, respectively, open the commutators, 
and perform the time integral. 



(Al) 



spin, and lead index. We introducing dot states |m) with 
eigenenergies and open the commutators. We only 
give the first of eight terms, the others are analogous: 



pt* = -I [Em - i?„) pTn - I dr trieads E ^ E E 

^ pq Qkcr i^u' 

(A3) 



Pmn ~ i [Em ^n) Pmn 

- ttE E ((z|(m|i/hybb)|j» ((j|(p|i/hyb|<7)K))p^°* KEp + - E, - e,) 

ij pq 

- W, ((z|(77l|i/hybb)|j))p^°*((j|(g|i^hyb|n)N)) 6iEq + - K - £,) 

- ((z|(m|i/hybb)|j))p^°*((j|(g|i^hyb|n)|z)) 5{Em + - Ep - e,) 
+ p;J,°*((z|(p|i/hyb|g)b-))((j|(g|i?hyb|n)|*)) 5{E„, + -E,- 



V PQ 



+*EE W» ((zi(m|i/hybb)ij)) mp\H^y^\qm)pt:' P 



1 



Ep + tj — Eq — ei 



W, ((z|(,7i|i/hybb)|j»Ppr(01(9l^hybl«)N» p 

W, ((i|(m|i/hybb)b-)>Pp°'(01('?l^hybb)N» P 



pr;((*|(p|i?hybk)|j))((j|('7|i?hybb)N» P- 



1 



■ Ei — Eq — €j 
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The second-order term contains the expression 



Tricads"ctk(T^ "ako-e Pleads 



(A4) 

All terms contain lesser or greater Green functions, G< 
or G^, respectively, which describe the lead correlations 
discussed in Sec. |lll Their Fourier transforms can be ex- 
pressed in terms of the Fermi function n p and the spec- 
tral function of the leads, AQ,k(j(w), as 



(A5) 



where /io, is the chemical potential of lead a. Performing 
the integral over r, we obtain 

"'1 akcr i^iy' 



npiuj - Ha) 



-LU + Ep-Eq- iO+ 



X {m\c^^\p){p\cl,Jq)p\ 



dot 
qn 



(A7) 



If we assume that the tunneling amplitudes tak^v = taau 
do not depend on the wave vector k, the sum over k 
can be performed, noting that the spin-resolved density 
of states is given by Daai^^) = {'2'^V)~^ ^^^Aakai^)- 
Here, V is the system volume. This leads to 



-dot i { W \ n 

Pmn — ~ ^n) Pmn 

- 71" ^ ^ ^ \^Daa{Ep - Eq- Pa) npiEp - Eg - p^) taauta^i,'{rn\Ci,a\p){j)\cl,„\q) + . . . | 

pq ocauv' 

-*^E E ^['^^-'^^'^"^''^A T-T+E WL,-Hc..|p)(p|ct,Jg)p^°* + ...|, (A8) 



where P indicates a principal value integral. 

We find that if the large-reservoir approximation is 
valid and if the tunneling amplitudes taka-u do not de- 
pend on k, only the lead density of states enters-^ Strong 
correlations in the leads, resulting in broad features in 
the spectral function, do no affect the tunneling. Under 
these conditions all information on k is lost so that it 
is only important whether states at a given energy (and 
a, a) exist. As discussed in Sec. IIIB[ broad features in 
the spectral function do not invalidate the Markov ap- 
proximation made here, since they correspond to rapid 
processes. 

APPENDIX B: FIRST-ORDER TERMS 

We briefly discuss first-order terms in the ME. Equa- 
tion (|65p shows that the only first-order term comes from 
the initial condition for Qpit^). The origin of the vanish- 
ing of first-order terms is that they contain equilibrium 



averages of single lead-fermion operators, cf. Eqs. ^ and 
(El. 

The absence of first-order terms might be surprising. 
Let us consider pdot with component p^* ^ 0, corre- 
sponding to the superposition of two states with electron 
number differing by one. If {n\H\-iy\^\m) 7^ 0, a single 
power of i/hyb is sufficient to lead to a change in p^l^m 
and p'f°l , analogous to the interaction of a superposition 
having an oscillating dipole moment with the light field. 
Superselection rule o^^i^^'^'' suggest that superpositions of 
states with different charge dephase so rapidly that they 
are unobservable. However, this does not help us at a 
fundamental level, since we would obtain the same equa- 
tion if the tunneling fermions were neutral. 

A first-order term is present if we expand the equa- 
tion for the full density operator, p{t) — e^*-^(*^*o) p[to), 
in powers of iZhyb- Since first-order processes thus exist, 
while Vp does not describe them, the required informa- 
tion must be contained in Qp(t). If they happen to be 
relevant, we must choose Qpito) 7^ 0. 
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